Electronegativity in quantum electronic transport 
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Electronegativity is shown to control charge transfer, energy level alignments, and electron cur- 
rents in single molecule tunnel junctions, all of which are governed by correlations contained within 
the density matrix. This is demonstrated by the fact that currents calculated from the one-electron 
reduced density matrix to second order in electron correlation are identical to the currents obtained 
from the Green's function corrected to second order in electron self-energy. 
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Prediction of electron transport across single molecules requires determination of electronic structure in the presence 
of open boundary conditions, whether using a non-equilibrium statistical or dynamical theory |l[ 0]. Statistical 
approaches concentrate directly on the non-equilbrium density matrix, whereas if the time evolution for a system 
driven from equilibrium is followed, attention is usually focused on the non-equilibrium Green's functions (NEGF) 
describing electron propagation. Treating electronic structure within transport theories requires an understanding of 
the intriguing and challenging to calculate effects of electron correlations. As exact approaches are limited to model 
systems or nanostructures with a small number of electrons, attention has focused on improving addition spectra both 
in the independent electron approximation 3( and by many-body treatments through the GW scheme [J]. In the 
following, we consider correlation corrections to independent particle models and relate conditions on the one-electron 
Green's function and reduced density matrix for calculation of currents within non-equilibrium theories. Correlation 
corrections to the density matrix are shown to correspond to improving ionization potentials (IPs) and electron 
affinities (EAs) given by Green's functions. This leads to a discussion of electron currents in terms of electronegativity: 
the impact of the electronegativity on charge transfer, energy level alignments, and current magnitudes is determined. 

Electron currents may be calculated from the one-electron reduced density matrix 5|i Ifil as 

^(r) = :^[Vr-Vr']p(r,r')|r'=r, (1) 
2.1 

with J the current density, r a position vector, and p the one-electron reduced density matrix (RDM); atomic units 
are implied unless otherwise given. As the current density operator is a one-body, non-local operator, it is clearly 
necessary to develop accurate approximations to the RDM to obtain reasonable results. From another viewpoint, 
calculation of the current can also proceed through computation of the one electron retarded and advanced Green's 
functions Gr,a and application of a Landauer-type formula [tJ: 

/= d^[/L(c^;ML)-/i?('^;Mi?)]Tr[rL(c^)Ga(a>)rfl,(^)AG,(cc>)], (2) 

with electron energy oj, r^^/j spectral densities, /l._r energy distributions with [i,h.R chemical potentials in the left (L) 
and right (R) electron reservoirs, and A is the correction due to correlations weighted by the spectral density of the 
electrodes and electron-electron spectral density on the molecule. The causal Green's function is related to the RDM 
via the relation 

p(r,r') = -^ </'dc.G(r,r';^), (3) 



with the complex integration performed along the Coulson contour. We begin by pointing out that the reduced density 
matrix obtained from a many-electron wavefunction corrected to second order in electron correlation is equivalent to 
the reduced density matrix arising from correcting IPs and EAs in the Green's function to second order in the electron 
self-energy ^. 

To proceed, the energy operator for a molecule within a tunnel junction is written in the form 

i?(A)= / dri/^t(i.)/i(r)^(r) + / dr dr'Vit(r)wHF(r, r')?A(r') 



-t-A[i J drdrV^(r)7/.t(r')w(r,r')?/'(r')i/'(r) dr dr'V>t (r)uHF(r, r')Vi(r')] , (4) 

with V the electron-electron interaction on the molecular region, whf the Hartree-Fock potential and 'tp^ip S'^e second 
quantized electron field operators. It is assumed that the Fock equations have been solved with electrode self-energies 
Y]LM. describe the interaction between the molecular region electrons and electrons in the reservoirs; external 
potentials are also included in the Fock operator. For A = 0, the Hamiltonian reduces to the Fock operator 

p 

with a\a creation and annihilation operators for Hartree-Fock states. For A = 1, the many-electron Hamiltonian is 
restored. A perturbation expansion in A is written for the many-electron wavefunction: 

l^r l^r(o) > +A|*(i) > -^A^l^-^^) > -H... . (6) 

For our choice of 0*'' order approximation, Brillouin's theorem insures that the first order wavefunction consists of only 
double electron excitations, on the other hand the second order term includes single through quadruple excitations. 
From 

p(r,r') =< ^'I^A^lrOV'WI* >, (7) 

. . n 

to first order in A the correction to the order density matrix vanishes [9]. The density matrix to second order is 

p«pW+AV'^- (8) 
The RDM may be represented as an infinite expansion over single electron states (j) 

p(r,r') = 5]pp,0;(r')0p(r). (9) 

pq 



Explicit calculation of the density matrix coefficients from eq. |6] through second order in A yields 



S l^y <abm><jk\\ab> ^ 



abk 



(ci + eb)(ej + <^k - fa - Cfc) 



1 < ij\\ac >< bc\\ij > 



< ij\\ac >< bc\\ij > 
. . {ei + Ej -ec- Ca) (ci + ej - Cc - £&) 

C 

1 < a6||ij >< aj||a6 > 1 < >< a5||zj > 



" 2^" (e, -e„)(e, + e, -e, -eb) 2 ^ (e, - e„)(e, + e, - e,, - e^) ^ ^ 

with < p(7||rs >=< pg|w|rs > — < pg|w|sr >. We use the convention whereby indices i,j,k,... label occupied, 
a,b,c, . . . label unoccupied, and p,q,r, . . . are used to label general (occupied or unoccupied) states in \'^^^'> >. 

Transmission resonances are given through the poles of the Green's functions and can be identified as IPs and 
EAs. Hence, it is reasonable to assume that if an independent particle picture is chosen to optimize IPs and EAs, it 
follows that prediction of currents from the NEGF approach will be improved. In this context, a model for transport 
is measured in terms of reproducing the molecular electronegativity. It is known that introduction of correlation 
corrections beyond independent particle models for the Green's function improves the predicition of IPs and EAs. 



The Green's function with second order self-energies has been studied by Pickup and Goscinski |g| leading to the 
following approximation 

1 ^ < ab\\pi >< qi\\ab > 1 ^ < ij\\pa >< qa\\ij > 

= [uJ - epjdpq - - y ■ -> ■ . (13) 

2 Lu + Si - ea - 2 ^ Lu + ea ~ ei - ej 

tab ija 

The lowest order improvement to Koopmans' IPs and EAs are obtained from the poles of the diagonal elements 
of G(uj). It is found the self-energy corrects Koopmans' IP through terms describing orbital relaxation and pair 
correlations; a similar interpretation holds for corrections to the EAs [8]. Within this approximation, it is also possible 
to determine the density matrix directly from eq. [31 the resulting density matrix coincides exactly with the density 
matrix calculated from eq. [6] through O(A^). Hence calculating the density matrix through second order in electron 
correlation and correcting IPs and EAs with second order self-energies S*-^-* will lead to the same predictions for 
electron current. For moderate electron correlations, improving spectra for independent particle models or explicitly 
including correlations in the RDM are equivalent. 

Recently a criterion for selecting an independent particle model for quantum electronic transport was given as the 
set of single particle states yielding an approximate density matrix with maximal overlap to the exact RDM jlo| . The 
single electron states diagonalizing the RDM are natural orbitals (NOs) pj| and their eigenvalues pt are known as 



natural occupations. If one asks what is the best finite expansion approximation p to the exact RDM 



^P = j \p-p?dTdr', 



(14) 



it is found that including the first n natural orbitals with the largest occupancies for a truncated expansion eq. [5] 
fulfills the least squares condition 9| . We consider the couplings between density matrix coefficients in eq. [10] by 



(15) 



writing 

Pi j Pi a 
Pai Pab 

with (ij), (a 6), and (ia) denoting occupied-occupied, unoccupied- unoccupied, and occupied- unoccupied spaces re- 
spectively, with occupations referred to the 0*'* order wavcfunction. The natural orbitals to second order in electron 
correlation are given by the eigenfunctions of eq. 1151 Constructing the "best" independent particle picture in the 
sense of eq. [14] implies occupying a single Slater determinant by the first n,, natural orbitals. We have previously 
shown numerically that a single determinant composed of the largest occupation number NOs can lead to essentially 
the same results as a full many-body treatment for tunneling through alkanes 10|. For a single determinant approx- 
imation, the density matrix is idempotent = p which occurs since the first rif. occupations are equal to 1 with all 
others 0. Hence a measure for the quality of a single determinant approximation is how well the eigenvalues of eq. [15] 
approximate the idempotency condition. As the pi a couplings between the occupied and unoccupied spaces becomes 
stronger, the occupations of the 0*'* order states can become significantly less than unity. From many-body theory 
it is well understood what this condition implies: a single determinant or independent particle picture is no longer 
useful as a 0*'' order wavcfunction. For weak to moderate correlations, the Green's function approach can achieve 
improved IPs and EAs by a low order approximation to the electron self-energy. As natural occupancies in the 0*'' 
order wavcfunction become very much less than unity, a perturbation expansion about an independent particle picture 
loses meaning and even higher order corrections to |vI/(") > will not correct IPs and EAs on the molecular region. In 
a similar context, this is seen as the failing of the GW approximation for systems with multi-determinantal ground 



states 12] or in strongly correlated electron transport 



For strong electron correlations coupled-cluster theory offers 



13|, 



a convenient nonperturbative framework from which higher order approximations to the density matrix follow 
alternatively correlated one particle methods to infinite order can be chosen to yield correct IPs and EAs. 

The MuUiken electronegativity given as {IP + EA)/2 is a useful measure of charge transfer, and it is charge transfer 



that determines molecular level alignments relative to electron reservoir energies 



15|. Predicting level alignments 



IJ. In 



correctly for molecules bonded between electrodes is essential for accurate current-voltage characteristics 
the Hartree-Fock approximation, charge transfer is under-estimated as hybridization to virtual states is weak. In 
the local density (LDA) and generalized gradient (GGA) approximations to density functional theory (DFT), charge 



transfer is over-estimated [17| . These effects are demonstrated for the case of hexenedithiol bonded between two gold 
clusters in fig. [T] where the highest occupied-lowest unoccupied energy gap in the molecular orbitals (HOMO-LUMO 
gap) and molecular electronegativity is given against charge transfer relative to molecular hexenedithiol. For a large 
HOMO-LUMO gap or weak electronegativity, charge transfer is small. For small HOMO-LUMO gaps typical of GGA 
and LDA, over-estimation of charge transfer is confirmed. Hybrid functionals correct the charge transfer to some 



extent, but this correction is not systematic 



17|. 



We introduce a simple correlated model for a molecular chain and investigate the effect of over- and under-estimation 
of electronegativity on electron transport. We use the following model Hamilitonian for an infinite chain: 

iJ= - 7L {clc^_j^+ h.c.) + {eL + VL)cic^-jLM{cLj)_3 + h.c.) 

n<— 3 n<— 3 

+3 

+ J2 (^^-f + ^" ^" " 7A/(&-3 b-2 + bli h + blh + h.c.) - Tm{V-2 ^-i + b\ h + h.c.) 

n— — 3 

- 7A/fl,(^l3C+4 + h.c.) + {eR + Vb) c\ c„ - 7^ ^ {c\ c^^^ + h.c.) (16) 

ri>+3 n>+3 

Six central sites of the chain are labelled -3, -2, -1, 1, 2 ,3 (i.e. there is no site) and are treated as the molecular region 
with 6^,6 creation and annihilation operators for electrons on the molecule. The electron reservoirs are described 
by the atomic sites extending towards the left and right away from the central molecular sites with creation and 
annihilation operators c^, c for the reservoir electrons. The site energies are given by — and eM for the reservoir 
and molecular regions, respectively. The volage applied across the molecular junction is described by the voltages 
Vl 7^ Vr in the reservoirs and the voltage drop Vn across the molecular sites is scaled linearly between the values Vl 
and Vr. The nearest neighbor interactions are = 1r within the electrode regions, and there are two molecular 
site-site interaction rj\/ and 7m representing single and double bonds, respectively, on the molecular region as a 
simple model for hexenedithiol, and 7la/ — Jmr. determine the molecule-electrode couplings. The eigenstates of 
the molecular Hamiltonian are found with the electron-electron self-energies and exact electrode self-energies are 
introduced describing coupling to the electrodes |l8| . The resulting single electron states are taken as the expansion 
functions for the correlated version of the model obtained from Hq Hq + v, with v the pairwise perturbation 
interactions about the mean field solution as in eq. [4]with A = 1. Current- voltage characteristics are calculated using 
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eq. [21 We use a simplified form of the self-energy such that the interaction matrix elements in eq. [13] are approximated 
as < mII^s >Ri U. In fig. [21 the HOMO-LUMO gap for the molecular region is given as a function of U demonstrating 
that the electronegativity on the molecular region may be systematically controlled through the electron-electron self- 
energies. The results for the current voltage characteristics from the model are presented in fig. [3l The independent 
particle or uncorrelated model occurs for U = and increasing U corresponds to increasing electron correlations on 
the molecular region. At U — 0, currents at low voltages are much lower than when the E*^^^ term is allowed to correct 
IPs and EAs; in this case, the highest lying occupied states are too low (IPs too high) and the lowest lying unoccupied 
single electron states are too high (EAs too low) with respect to the Fermi level. Under these conditions neither 
occupied or unoccupied states enter into the voltage bias window at low voltages, and this level of electronic structure 
treatment corresponds to a Hartree-Fock approximation. Increasing correlations on the molecular region, the highest 
occupied states near the Fermi level enter the bias window at lower values of voltage, followed by the introduction of 
the unoccupied states at higher voltage bias (this sequence is due to the relative position of the Fermi level relative 
to occupied and unoccupied states for this example). The correlations on the molecular region serve to shift up 
occupied levels relative to the Fermi level leading to reduced IPs, whereas increasing correlations systematically lower 
the lowest lying unoccupied states leading to increased EAs. Increasing correlations continue to reduce the IPs and 
increase EAs until eventually electronegativity is under-estimated. The impact on the current-voltage characteristics 
is that the molecular levels enter the bias window at very low values of applied voltage resulting in large current 
magnitudes. Larger values of U correspond to the use of LDA or GGA exchange-correlation potentials within DFT 



where the strong over-estimation of charge transfer is known to occur 17|. Hartree-Fock and Kohn-Sham (using 
LDA or GGA) orbitals are not appropriate independent particle models for electron transport due to strong under- 
and over-estimation of charge transfer, respectively. The results of fig. [H clearly show the impact on current voltage 
characteristics for these two extremes. 

Correcting electronegativity is equivalent to maximizing overlap to the reduced density matrix: this is true to 
low orders in electron correlation and of course the correct electronegativity and density matrix are found at the 
exact many-body solution. In general, improving descriptions for the RDM and electronegativity with the methods 
described will lead to improved prediction of electron currents in systems with moderate electron correlations. The 
best independent particle picture within this context is a single determinant comprised of natural orbitals; any attempt 
to refine single electron models for transport should lead to electron wavefunctions that approximate natural orbitals. 
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In the case of Green's function approaches, moderate electron correlations imply the need to include electron-electron 
self-energies to describe quasi-particle propagation. For strong correlations, a single determinant wave function is not 
an adequate approximation to predict IPs and EAs and perturbation corrections about a single reference state fail- 
thus complicating treatment of molecular junctions with Green's function approaches. However, in all cases, from 
weak to strong correlations, the criterion to maximize overlap to the exact reduced density matrix leads to improved 
predictions for electron currents. 
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FIG. 1: Highest occupied and lowest unoccupied (HOMO-LUMO) energy gap and electronegativity versus charge transfer for 
hexenedithiol bonded to linear gold chains. Calculations have been performed with the TURBOMOLE program system [19 1. 
All calculations have been performed using the auc-cc-pVDZ basis set for carbon 2^ and split valence polarized valence basis 
for all other atoms, including a sixty electron effective core potential for the gold atoms Calculations have been performed 
using the Hartree-Fock and density functional theory calculations using hybrid (B3-LYP), generalized gradient approximation 
(GGA/PBE), and local density approximation (LDA/PW) exchange-correlation functionals. 
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FIG. 2: HOMO-LUMO gap for the model system defined by eq. [16] as a function of the electron-electron self energy as varied 
through the interaction parameter U . The reduction in the gap demonstrates the effect of electron-electron self-energy on the 
molecular electronegativity. 
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FIG. 3: Current voltage characteristics for the model Hamiltonian of eg. 1161 em ~ l.OeV, eL ~ ~ epormi = 0.0, 7a/ = 4.54eV^, 
Tm = l.beV, 7l = 7fl = lO.OeV^, 7lm = Jmr = 2AeV . Electronegativity is modified by varying U , with values as labeled 
within the figure. 



